Draft version December 6, 2010 

Preprint typeset using LAT^X style cmulatcapj v. 11/10/09 



THE BLACK HOLE MASS IN BRIGHTEST CLUSTER GALAXY NGC 6086 

Nicholas J. McConnell 1 , Chung-Pei Ma 1 , James R. Graham 1,2 , Karl Gebhardt 3 , Tod R. Lauer 4 , Shelley A. 

Wright 1 , and Douglas O. Richstone 5 

Draft version December 6, 2010 

ABSTRACT 

We present the first direct measurement of the central black hole mass, Mm, in NGC 6086, the 
Brightest Cluster Galaxy (BCG) in Abell 2162. Our investigation demonstrates for the first time that 
stellar dynamical measurements of M, in BCGs are possible beyond the nearest few galaxy clusters. 
We observed NGC 6086 with laser guide star adaptive optics and the integral-field spectrograph 
(IFS) OSIRIS at the W.M. Keck Observatory, and with the seeing-limited IFS GMOS-N at Gemini 
Observatory North. We combined the IFS data sets with existing major-axis kinematics, and used 
axisymmetric stellar orbit models to determine M, and the i?-band stellar mass-to-light ratio, 
M*/L fl . We find M. = 3.6^ x 10 9 M Q and M+/L R = 4.6l°;? M /L© (68% confidence), from 
models using the most massive dark matter halo allowed within the gravitational potential of the 
host cluster. Models fitting only IFS data confirm M. ~ 3 x 10 9 M Q and M*/L fl ~ 4 M o /L , with 
weak dependence on the assumed dark matter halo structure. When data out to 19 kpc are included, 
the unrealistic omission of dark matter causes the best-fit black hole mass to decrease dramatically, 
to 0.6 x 10 9 M Q , and the best-fit stellar mass-to-light ratio to increase to 6.7 Mq/Lq^. The latter 
value is at further odds with stellar population studies favoring M*/L^ ~ 2 M Q /L Q . Biases from 
dark matter omission could extend to dynamical models of other galaxies with stellar cores, and 
revised measurements of M, could steepen the empirical scaling relationships between black holes 
and their host galaxies. 



1. INTRODUCTION 

It is increasingly accepted, both observationally 
and theoretically, that supermassive black holes 
are ubiquitous at t he centers of elliptical galax- 
ies (jMagorrianl I1998D . The black hole mass, M., 
correlates with various host properties, including 
bulge luminosity, L (e .g., iKormendv fc R ichstone fl995t 
IMarconi fc Hun t 2003), and stellar velocity dispersion, a 
(e.g.. IFerrarrese fe Merrittl [2000l : IGebhardt etaH [2000 1. 
These empirical correlations have been established from 
approximately 50 galaxies in which M, has been de- 
termined from motions of stars, gas, or masers under 
the direct gravitational influence of the central black 
hole. Although galaxies with Ly ~ 10 9 — 10 11 Lm y are 
well-represented in th is sample (e.g.. lHaring fc Rbal2004l : 
iGiiltekin et "all 120091 ) . there are very few measurements 
of M, in the most luminous galaxies. 

Brightest Cluster Galaxies (BCGs) are the most lumi- 
nous galaxies in the present-day universe (Ly ~ 10 10 5 — 
10 115 L Qi y). Direct measurements of M» in these galax- 
ies have been lacking because very few kinematic stud- 
ies spatially resolve the black hole radius of influence, 
»"mf = GM./ct 2 . The M, — a relation predicts typical val- 
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ues of T"i n f ~ 30 pc in BCGs; predictions from the M, — L 
relation are a few times larger. BCGs' low central sur- 
face brightnesses exacerbate the challenge of obtaining 
high-quality stellar absorption spectra at angular scales 
comparable to n n f . To date, stellar dynamical measure- 
ments of M, in BCGs have been limited to the nearest 
groups and clusters: M87 in Vi rgo (e.g., iSargent et al.l 
119781 . IGebhardt fc Thomas! [20091) and NGC 1399 in For- 
nax (|Houghton et al J l2006t IGebhardt et al J 120071 ). In a 
small number of BCGs, M, can be measured from emis- 
sion line kinematics in a resolved disk of ionized gas. 
IDalla Bonta et al.l ([2009D have used STIS on the Hubble 
Space Telescope (HST) to examine disks at the centers 
of three BCGs beyond 50 Mpc, reporting two measure- 
ments of M m and one upper limit. 

BCGs are distinct from other giant elliptical galaxies 
in several respects. Two such distinctions are partic- 
ularly intriguing with regards to the evolutionary con- 
nections between galaxies and their central black holes. 
First, BCGs are preferentially found near the gravita- 
tional centers of galaxy clusters, where cosmological dark 
matter filaments intersect. Second, BCG luminosities 
vary more steeply with a than the canonical L oc cr 4 rela - 
tionship for elliptical galaxies (|Oegerle fc Hoe sscl 1991); 
lLauer et al.l ([2007; have found L oc a 7 for BCGs and 
other core-profile galaxies. The steep relationship be- 
tween L and a in very massive galaxies requires one or 
both of the M, — a and Af, — L relationships to diffcr 
from the correlations observed in lower-mass galaxies. 
Direct measurements of M, in a statistically significant 
sample of BCGs will reveal the forms of these relation- 
ships for the most massive galaxies, and will help dis- 
criminate different evolutionary scenarios fo r BCG s. For 
instance, iBovlan-Kolchin. Ma fc Quataertl (|2006l ) have 



demonstrated that gas-poor galaxies merging on radial 
orbits could produce the steep relation between L and 
a. With little gas available for star formation or black 
hole accretion, the remnant galaxy and black hole would 
remain on the same Af, — L relation as the progeni- 
tors. These radial mergers could occur at the intersec- 
tion of cosmological fil a ments . In one counterexample, 
iRuszkowski fe Springell ()2009f ) performed a zoom-in res- 
imulation of a single 10 15 M Q galaxy cluster selected from 
a cosmological iV-body simulation, and produced a BCG 
that remained on nearly the same L — a relation as the 
fainter galaxies. A larger sample of resimulated clusters 
would help assess the relative frequency of radial orbits 
and their impact on the scaling relations of BCGs. Alter- 
native scenarios f or BCG growth, such as early- time ma- 
jor mergers (e.g.. iMerrit t 1985; Trcmainc 199( J) or "can- 
nibal i sm" of smaller galaxies (e.g ., lOstriker fc Tremaind 
119751 : lOstriker fc Hausmanll 19771) . potentially could pro- 
duce lower values of AT., matching predictions from the 
M. — a relationship. In these scenarios, the final black 
hole mass could depend upon a number of factors, such 
as the orbits, gas fractions and disk-to-bulge ratios of 
merging galaxies. 

In addition to providing clues toward BCG evolu- 
tion, empirically establishing the high-mass forms of the 
M. — g and M m — L relationships will provide new con- 
straints for the number density of the Universe's most 
massive black holes. The most luminous high-redshift 
quasars a re inferred to host black holes exceeding 10 10 
Mfp ( e.g., lBechtold et al.ll2003t lNetier]|2003t FVestergaardl 
120041 ) . but thus far no such objects have been detected 
in the local Universe. BCGs in nearby Abell clusters 
potentially c ould host black holes with Af, > 10 9 ' 5 M Q 
(|Lauer et aI1l2007D . 

Another motivation for measuring Af, in BCGs is that 
the faint centers of these galaxies likely arise from "core- 
scouring," whereby stars are ejected from the galatic cen- 
ters by an in-spiraling pa ir of supermassive black holes af- 
ter a major merger (e.g.. lEbisuza ki. Makino fc Okumural 
I1991D . Given theoretical expectations for the efficiency of 
core scouring, a galaxy's past merger history can be esti- 
mated by comparing M, to the total luminosity deficit i n 
the core (|Lauer et al.ll2007t iKormendv fc Benderll2009f ). 

In this paper, we report measurements of Af. and the 
f?-band stellar mass-to-light ratio, M^/L^, in NGC 6086, 
the BCG of Abell cluster 2162. Our investigation marks 
the first attempt to measure Af, using stellar dynam- 
ics in a BCG beyond Virgo. Future papers will describe 
measurements of Af, in a larger sample of BCGs. For 
BCGs at ~ 100 Mpc, 8- to 10-meter telescopes with 
adaptive optics (AO) are required to obtain good spectra 
on ~ 0.1" spatial scales. Laser guide star adaptive optics 
(LGS-AO) enables the study of targets without a bright 
nearby guide star. We use integral-field spectrographs 
(IFS) to obtain full 2-dimensional spatial coverage, which 
places tighter constraints on stellar orbits. Our orbit 
models include a dark matter co mponent in the gravi- 
tation al potential, as described in iGebhardt fc Thomas! 
(2009). In this paper, we emphasize methods for pair- 
ing IFS data with axisymmetric orbit models, and for 
assessing errors in Af, and M^/L^. 

NGC 6086 is a cD galaxy at the center of Abell 
2162. Like many BCGs, it exhibits radio emission 



(|Ledlow fc Owenl I1995D , likely from low-level accretion 
onto the central black hole. We have derived an effec- 
tive stellar velocity dispersion of 318 km s _1 in NGC 
6086, using meas urement s from Carter, Bridges & Hau 
(1999; hereafter ICBH99D . This would correspond to 
a black hole mass of 9 x 10 s M Q , if NGC 6086 were 
to fol low the mean Af, — a relation of Giiltekin et al. 
(p00l . The F-band luminosity of NGC 6086 is 1.4 x 10 11 
L 0j y, from M v = -23.11 in lLauer et al.l (|2007h : the 
corresponding bla ck hole mass predicte d from the mean 
M.-L relation of lGultekin et al.l (|2009l ) is 1.3 x 10 9 M . 
Abell 2162 is a relatively small galaxy cluster at redshift 
z = 0. 032, with a richness class of based on 37 mem- 
bers (jAbell. Corwin fc QlowinlH989[ ). and a line-of-sight 
velocity dispersion of 302lgg 2 km s -1 (Zabludof f et al.1 
[1991 . NGC 6086 is offset from the average radial veloc- 
ity of Abell 2162, by 82 km s" 1 (jLaine et al.l 120031 ). 

This paper is organized as follows. In §2, we describe 
our photometric data of NGC 6086, our IFS observations 
at Keck and Gemini observatories, and the subsequent 
data reduction procedures. In §3, we describe our proce- 
dures for extracting 2-dimensional kinematics from IFS 
data, and compare our resulting measurements in NGC 
6086 with other studies. In §4 we review the stellar orbit 
modeling procedure. We also report our measurements of 
Af. and M^/L^, and describe how these measurements 
depend on the assumed dark matter halo profile. We 
estimate confidence intervals for Af, and M^/L^, and 
discuss both tested and un-tested systematic errors. In 
§5 we compare our results to predictions from the Af, — a 
and Af, — L relationships, and discuss whether the effect 
of dark matter on stellar orbit models of NGC 6086 can 
be generalized to reveal biases in measurements of Af, in 
other galaxies. Appendix A contains a detailed descrip- 
tion of systematic errors from stellar template mismatch 
and uncertain PSFs. Appendix B contains our full set of 
measured line-of-sight velocity distributions (LOSVDs). 

Throughout this paper, we assume ffo = 70 km s , 
TO = 0.27, il\ = 0.73, and an angular-diameter distance 
of 133 Mpc to NGC 6086. One arc sec corresponds to 
0.64 kpc at this distance; for a — 318 km s , n n f = 
0.066" x (Af./10 9 M Q ). 

2. OBSERVATIONS 
2.1. Photometry 

We use a combination of f?-band (0.6 /jm) and f-band 
(0.8 /im) photometry to constrain the stellar mass pro- 
file of NGC 6086. For radii out to 10" we adopt the 
high-resolut i on su rface brightness profile presented in 
lLaine et al.l (|2003l ). obtained with WFPC2 on the Hubble 
Space Telescope. This surface brightness profile has been 
corrected for the WFPC2 point-spread function (PSF) 
by applying the Lucy-Richard son deconvolution method 
()Richardsonlll972tlLucvlll974t ) ; specific details of the im- 
plementation are described in lLaine et al.l (|2003t ). 

At larger radii out to 86" we use f?-band data from 
Lauer, Postman & Strauss (private communication), ob- 
tained with the 2.1-m telescope at Kitt Peak National 
Observatory (KPNO). The KPNO data have a field-of- 
view (FOV) of 5.2' x 5.2', which enables accurate sky 
subtraction. To create a single surface brightness pro- 
file, we assessed the individual profiles from WFPC2 and 
KPNO data at overlapping radii between 5" and 10". 



2 



10 2 
10° 



- 10" 



10' 
10 s 





1 ' ' 

""^^"^ WFPC2 


1 

T 
T 

1 
1 






1 

KPNO 
T 
1 

V 1 

\ 1 


1 < < 


1 I I 


i -X 
i **.\ 
i "%\ 
i **,\ 
i \\ 
i i i i-\ i 



0.2 



0.4 



Flux Relative to Peak 
0.6 



0.8 



1.0 



0.1 



1.0 10.0 

r (arc sec) 



100.0 



Fig. 1. — De-projected ij-band stellar luminosity density vs. ra- 
dius along the major axis (solid line) and minor axis (dotted line) of 
NGC 6086. The dashed vertical lines mark the outermost extents 
of photometric data from HST /WFPC2 and KPNO. Luminosity 
densities beyond 86" are derived from a de Vaucouleurs surface 
brightness profile. The H-band surface brightness at 86" is 24.3 
mag arcsec - 2 on the major axis. 

We measured the average R — I color for these radii 
and added it to the WFPC2 profile. The two profiles 
were then stitched together such that their respective 
weights varied linearly with radius between 5" and 10": 
the WFPC2 data contribute 100% to the combined pro- 
file for r < 5" and the KPNO data contribute 100% for 
r > 10". PSF deconvolution was not necessary for the 
KPNO data, as they contribute to the combined surface 
brightness profile at radii well beyond the seeing full- 
width at half- maximum (FWHM). Our translation of the 
WFP C2 profile to J?-ban d assumes no R — I color gra- 
dient; lLauer et al.l (|2005l ) find a median color gradient, 
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of -0.03 magnitudes for BCGs and other core 
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profile galaxies. 

At radii beyond 1", isophotes of NGC 6086 all have 
major-axis position angles within 5° of true north, with 
an average apparent axis ratio of 0.7. We adopt 0° east 
of north as the major-axis position angle of NGC 6086, 
and we assume edge-on inclination. We deprojected the 
surfac e brightness using the procedure of lGebhardt et al.l 
(|1996f) . which assumes spheroidal isodensity contours. 
The resulting major- and minor-axis luminosity density 
profiles are shown in Figure [TJ 

2.2. Spectroscopy 

We performed integral-fi eld spectroscopic ob servations 
of NGC 6086 with OSIRIS (jLarkin et alJl20 06) on the 10- 
m W. M. Keck II telescope and GMOS-North (Allington- 
Smith et al. 2002; Hook et al. 2004) on the 8-m Gemini 
Observatory North telescope. The instrument properties 
and our observations are summarized in Table [T] Our 
observations with OSIRIS used the W. M. Keck Obser- 
vatory laser guide star adaptive optics (LGS-AO) system 
(Wizinowich et al. 2006; van Dam et al. 2006); the inner 
component of the resulting .ff-band (1.6 /im) PSF has an 
FWHM value of w 0.1". The CMOS data were collected 
under excellent seeing conditions; images of point sources 
from the Gemini North Acquisition Camer£[j indicate an 

6 http:/ /www. gemini.edu/sciops/telescopes-and- 

sites / acquisition- hardware- and-techniques / acquisition-cameras 




Fig. 2.— Total H-band flux for NGC 6086, using collapsed spec- 
tra from OSIRIS. The dashed line traces the major axis of the 
galaxy, with the arrow pointing north. Thick black lines enclose 
the spatial bins corresponding to the spectra displayed in Figure \3\ 

J-band FWHM of 0.4". 

In Figure [2] we display the reduced mosaic of NGC 
6086 from OSIRIS, summed over all spectral channels. 
Usable data from OSIRIS and GMOS extend to radii 
of 0.84" and 4.9", respectively. For radi i out to 30" we 
use major-axis kinematics from ICBH99L obtained with 
the Intermediate Dispersion Spectrograph (IDS on the 
2.5-m Isaac Newton Telescope. 

2.2.1. OSIRIS 

OSIRIS is a near- infrared (NIR), integral- field spec- 
trograph built for use with the Keck AO system. It fea- 
tures 2-dimensional spatial sampling at four scales be- 
tween 0.02" and 0.1". We observed NGC 6086 with the 
0.05" spatial scale, which provided adequate signal-to- 
noise and placed several pixels within the radius of influ- 
ence. To minimize noise in individual spectra, we used 
the broad ii-band filter, which covered several Ais = 3, 
12 CO bandheads at observed wavelengths from 1.54 /im 
to 1.71 /im (at z « 0.032). We chose to detect ff-band 
features instead of the more prominent v =0-2 12 CO 
bandhead in K — band, which suffered from higher ther- 
mal background at an observed wavelength of 2.37 /im. 

We recorded 9 science exposures of the galaxy center 
and 5 sky exposures of a blank field 50" away, for to- 
tal integration times of 2.25 hr and 1.25 hr, respectively. 
Our dithers repeated an "object-sky-object," sequence, 
such that every science frame was immediately preceded 
or followed by a sky frame. We also recorded spectra of 
9 spectral template stars, using the same filter and spa- 
tial scale as for NGC 6086. To measure telluric absorp- 
tion, we recorded spectra of several A0V stars, covering 
a range of airmasses similar to those for NGC 6086 and 
template stars. 

We used version 2.2 of the OSIRIS data reduction 
pipeline^ to subtract sky frames, correct detector arti- 
facts, perform spatial flat-fielding, calibrate wavelengths, 
generate data cubes with two spatial dimensions (x,y) 
and one spectral dimension (A), and construct a mosaic 
of NGC 6086 from multiple data cubes. The pipeline 
uses an archived calibration file to perform spectral ex- 



7 http://www.ing.iac.es/Astronomy/instrumcnts/ids/ 

8 available from the UCLA Infrared Laboratory, 
http: / /irlab. astro. ucla.edu/osiris/pipeline. html 
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TABLE 1 

Summary of spectroscopic observations 



Instrument 


UT Date 


A Range 


AA 


Ax 


tint 


PA 


FWHM 






(nm) 


(nm) 


1 arr cpf 1 

V £11 DCL 1 


w 






(1) 


(2) 


(3) 


(4) 


(5) 


(6) 


(7) 


(8) 


OSIRIS 


May 13-14, 2008 


1473 - 1803 


0.2 


0.05 


9 x 900 


-45 


0.10 


GMOS 


April 25, 2003 


744 - 948 


0.1377 


0.2 


5 x 1200 





0.4 


IDS 


June 10-16, 1996 


493 - 573 


0.2 


0.4 x 3.0 


22,800 





1.3 



Notes: Column 1: instrument. OSIRIS (OH-Suppressing Infra-Red Imaging Spectrograph) was used on Keck II with LGS-AO. 
GMOS (Gemini Multi-Object Spectrograph) was used on Gemini North. IDS (Inte rmediat e Dispersion Spectrograph) was used 
on the Isaac Newton Telescope; here we summarize the observations published by C BH 99L Column 2: date(s) of observations. 
Column 3: observed wavelength range. Column 4: spectral pixel scale in 3-d data cubes, for OSIRIS and GMOS data. FWHM 
spectral resolution, for IDS data. Column 5: angular spacing of lenslets, for OSIRIS and GMOS data. Pixel scale along slit 
x slit width, for IDS data. Column 6: number of science exposures x integration time per exposure. For IDS data, the total 
integration time of 6.33 hours is reported from CBH99. Column 7: position angle of the long axis for OSIRIS and GMOS, or 
the slit for IDS, in degrees east of north. Column 8: PSF FWHM at science wavelengths. 
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Fig. 3.— Top: OSIRIS spectrum of the center of NGC 6086 
(0.01 arcsec 2 ; S/N = 21). Middle: OSIRIS spectrum 0.49" from 
the center of NGC 6086 (0.28 arcsec 2 ; S/N = 39). The dashed and 
dotted portions of OSIRIS spectra are excluded from the kinematic 
fitting, with red dotted portions specifically indicating regions of 
telluric OH contamination. Each thick blue line is the M4III tem- 
plate spectrum from OSIRIS, convolved with the best-fit LOSVD 
for the respective galaxy spectrum. Bottom: Spectrum of template 
star HD 110964 (M4III), from observations with OSIRIS. This is 
the only star used in our final extraction of kinematics from OSIRIS 
data. 

traction of the raw spectra across the detector and as- 
semble a data cube; the calibration file was generated by 
illuminating individual columns of the OSIRIS lenslet ar- 
ray with a white light source. We used custom routines 
to remove additional bad pixels from detector images, 
extract 1-dimensional stellar spectra from 3-dimensional 
data cubes, and calibrate galaxy and template spectra 
for telluric absorption. Although 1-dimensional stellar 
spectra from OSIRIS comprise an average over many spa- 
tial pixels, spatial variations in instrumental resolution 
are negligible relative to the velocity broadening in NGC 
6086: (Acr inst ) 2 - 5 x l(T 3 o- 2 . 

Contamination from telluric OH emission presents a se- 
vere challenge for observing faint, extended objects with 
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Fig. 4.— Top: GMOS spectrum of the center of NGC 6086 
(0.24 arcsec 2 ; S/N = 91). Dashed portions are excluded from the 
kinematic fitting. Middle: G9III template spectrum from GMOS, 
convolved with the best-fit LOSVD. Bottom: Spectrum of template 
star HD 73710 (G9III), from observations with GMOS. This is the 
only star used in our final extraction of kinematics from GMOS 
data. 

OSIRIS. The small field of view (0.8" x 3.2" for broad- 
band observations at 0.05" per spatial pixel) does not 
allow for in-field sky subtraction, and subtracting consec- 
utive science and sky frames only provides partial correc- 
tion, as the relative flux from different vibrational transi- 
tions in OH varies on timescales of a few minutes. After 
subtracting a sky frame from each science frame, we are 
forced to discard the spectral channels with strong resid- 
ual signals from OH, which compose approximately 15% 
of our spectral range. In Figure[3J we illustrate represen- 
tative spectra from OSIRIS and distinguish kinematic 
fitting regions from residual telluric features. At both 
ends of the ii-band spectrum, atmospheric water vapor 
acts as an additional contaminant. We have restricted 
our kinematic analysis to observed wavelengths between 
1.48 and 1.73 /urn. 

A second challenge for studying the centers of galax- 
ies with OSIRIS is accurate determination of the PSF. 
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We must construct an average PSF for mosaicked data 
from several hours of observations over multiple nights, 
during which seeing conditions and the quality of AO cor- 
rection can change significantly. To estimate the PSF, we 
recorded a one-time sequence of exposures of the LGS- 
AO tip/tilt star for NGC 6086, using the OSIRIS spec- 
trograph with the same filter and spatial scale settings as 
for our science frames. Data cubes were then collapsed 
along the spectral dimension to produce images of the 
star. In Appendix [3] we discuss different methods for 
estimating the PSF, and how PSF uncertainty influences 
our modeling results. 

2.2.2. GMOS 

GMOS-N is a multi-purpose spectrograph on Gemini 
North. GMOS includes an IFS mode, in which hexag- 
onal lenslets divide the focal plane and fibers map the 
2-dimensional held to a 1-dimensional slit configuration. 
A second set of lenslets samples a field ~ 60" from the 
science target, allowing for simultaneous sky subtraction. 
We observed the center of NGC 6086 in IFS mode with 
the detector's CaT filter, detecting the infra-red Call 
triplet at observed wavelengths between 0.87 and 0.90 
/im. A representative spectrum from the center of NGC 
6086 is shown in Figure |H GMOS data were reduced 
using version 1.4 of the Gemini IRAF software pack- 
agfl This standard pipeline subtracts bias and over- 
scan signals, removes cosmic rays, mosaics data from 
three CCDs, extracts spectra, corrects throughput vari- 
ations across fibers and within individual spectra, cali- 
brates wavelengths using arc lamp exposures, computes 
an average sky spectrum, and performs sky subtraction. 
We stored individual spectra from each GMOS exposure, 
along with their spatial positions relative to the center 
of NGC 6086, for eventual spatial binning. 

With seeing-limited spatial resolution, GMOS poorly 
resolves the black hole sphere of influence. Nonetheless, 
kinematics derived from GMOS provide a good com- 
plement to those from OSIRIS. The Call triplet region 
in GMOS spectra has a more clearly- defined continuum 
than -ff-band spectra, and with less telluric contamina- 
tion, as is evident from comparing Figures [3] and [4] As a 
result, line-of-sight velocity distributions (LOSVDs) ex- 
tracted from GMOS spectra have lower systematic errors 
than LOSVDs extracted from OSIRIS spectra. Addi- 
tionally, the angular region yielding high signal-to-noise 
spectra from GMOS is four times larger than that for 
OSIRIS. 

3. KINEMATICS 

Our dynamical models fit weighted and superposed 
stellar orbits to LOSVDs extracted from spectroscopic 
data. For both OSIRIS and GMOS data, we extract 
LOSVDs with a Maximum Penalized Likelihood (MPL) 
technique, which fits an LOSVD-convolved stellar tem- 
plate to each galaxy spectrum. The LOSVDs are non- 
parametric, defined at 15 radial velocity bins in our or- 
bit models. Representative LOSVDs from the central 
OSIRIS and GMOS bins are shown in Figure El and the 
full sets of LOSVDs extracted from OSIRIS and GMOS 
spectra are presented in Appendix |5] The MPL fitting 

9 available from Gemini Observatory, at 

http: / / www.gemini.edu/sciops / data-and-results / processing-software 
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Fig. 5.— Sample LOSVDs for NGC 6086. Left: LOSVD ex- 
tracted from the central spatial region, measured with OSIRIS 
(0.1" X 0.1"). Right: LOSVD extracted from the central spatial 
region measured with GMOS (~ 0.55" diameter). The solid blue 
line in each figure is the corresponding LOSVD generated by the 
best-fitting orbit model with the maximum-mass LOG dark mat- 
ter halo (M*/L fl = 4.7 M o /L ; M. = 3.5 X 10 9 M ; v c = 500 
kpc; r c = 8.0 kpc). The dotted red line in each figure is from the 
best-fitting orbit model with no dark matter halo (M*/Lfl = 6.8 
M /L ; Mm = 3.2 X 10 8 M ). For each LOSVD, Ax 2 is the dif- 
ference in the x 2 statistic for the two models: A% 2 > indicates 
that the model including dark matter (solid blue line) is a better 
fit. The full sets of LOSVDs from OSIRIS, GMOS, and ICBH99I 
are shown in Figures 1 171 through 1191 

method is described in det ail inlGebhardt et al.l ((2000b) , 
IPinknev et all (pOOl . and lNowak et al.l (|2008D T Here we 
describe the specific adjustments made for IFS data of 
NGC 6086. 

3.1. Extracting Kinematics from IFS Data 

In order to attain sufficient signal-to-noise (S/N) for 
effective kinematic extraction, we perform spatial bin- 
ning on our 2-dimensional grids of spectra from OSIRIS 
and GMOS. For the mean-normalized galaxy spectrum 
Y, mean- normalized stellar template T, and LOSVD C 
from the best fit over N c spectral channels, we define: 

S/N^ ( ^ [Yi — (T* £)J 2 / N c \ (1) 

At the very center of NGC 6086, we spatially bin spectra 
until S/N > 20 is achieved. This requires binning 2x2 
spatial pixels from OSIRIS; consequently our kinematic 
data have a central spatial resolution of 0.1", similar to 
the PSF FWHM. At the center of the GMOS mosaic, 
we bin 7 hexagonal pixels, corresponding to an approx- 
imate diameter of 0.55". The remaining spectra from 
each dataset are grouped to match angular and radial 
bins defined within the orbit models, and to maintain 
S/N between 25 and 40. Our resulting binning schemes 
for both OSIRIS and GMOS use only two angular bins 
on each of the positive (north) and negative (south) sides 
of the major axis. The angular bins span — 36.9° and 
36.9 — 90° from the major axis. Axisymmetric models 
perform LOSVD fitting in one quadrant of the projected 
galaxy. Symmetry about the major axis is enforced by 
co-adding spectra from the positive and negative (east 
and west) sides of the minor axis, before LOSVD ex- 
traction. LOSVDs extracted from the negative (south) 
side of the major axis are inverted before being input to 
the models. We define systemic velocity relative to the 
template star separately for OSIRIS and GMOS data. 
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Additionally, a spectral binning factor is necessary 
to smooth over channel-to-channel noise in spectra of 
NGC 6086. Our final kinematic extraction uses smooth- 
ing factors of 30 and 12 spectral pixels for OSIRIS and 
CMOS spectra, respectively. These values are chosen by 
comparing the best-fit LOSVDs from a large range of 
smoothing factors, and identifying the minimum factor 
above which LOSVDs in each dataset are stable between 
—500 km s _1 and 500 km s _1 . Our smoothing values 
are consis tent with the range of optimal values deter- 
mined by iNowak et al.l ()2008l ) for near-infrared spectra 
with S/N - 25 - 50. 

ii-band spectra from OSIRIS contain several absorp- 
tion features that are potentially useful for kinematic 
extraction, but some are compromised by incompletely 
subtracted telluric OH emission lines, which are masked 
from the fit. Three broad features are relatively insensi- 
tive to the narrow OH lines: the v = 3-6 12 CO bandhead 
at 1.6189 ,um rest, the v = 4-7 12 CO bandhead at 1.6401 
jum rest, and the Af = 2 band of OH between 1.537 
and 1.545 /im rest. Additionally, the Av = 2 OH band 
between 1.526 and 1.529 /mi rest does not intersect any 
strong telluric emission features. We have verified that 
these four spectral features offer a robust comparison be- 
tween stellar and galaxy spectra by repeating the fits 
with a large range of spectral smoothing factors. When 
we add other features to the fit, the root-mean-squared 
residual (RMS, essentially S/N^ 1 ) becomes unstable to 
small changes in spectral smoothing. 

To extract LOSVDs from CMOS spectra, we analyzed 
the A8498 and A8542 lines of Call. The third line in the 
well-known calcium triplet, A8662, is compromised by a 
flat-field artifact and discarded from kinematic analysis. 
The thick black lines in Figure[3]and Figure|4]indicate the 
spectral channels used in our final extraction of LOSVDs 
from OSIRIS and GMOS, respectively. 

Uncertainties for each LOSVD are determined by 100 
Monte-Carlo trials. In each trial, random noise is added 
to the galaxy spectrum, according to the RMS value of 
the original fit, and the LOSVD fitting process is re- 
peated. At each velocity bin, the uncertainties ac+ and 
<7£_ are computed from the distribution of trial LOSVD 
values. We then adjust the uncertainties in the wings of 
each LOSVD, so that C - a C - = 0. 

Stellar template mismatch can be a major s ource of 
systematic error in determining LOS VPs (e.g ., [CBH99; 
ISilge fc Gebhardtl l200l lEmsellem et~aT1 l200l . Toad- 
dress this issue, we have directly observed a diverse set 
of late-type template stars. Our nine templates from 
OSIRIS are giant, supergiant, and dwarf stars with spec- 
tral types from G8 to M4. We have found the most ap- 
propriate template for OSIRIS spectra of NGC 6086 to be 
HD 110964, an M4III star. In Appendix E we describe 
our method for choosing the optimal template star, and 
how fitting LOSVDs with different template stars influ- 
ences measurements of M, and M*/L#. We fit GMOS 
spectra of NGC 6086 with a single G9III template star, 
HD 73710, which we observed with GMOS. The calcium 
triplet region is less sensitive to template mismatch than 
other optic al and near-infrared regions used to measure 
kinematics (jBarth. Ho fc Sargentir2002l ). 

3.2. Two- Dimensional Kinematics in NGC 6086 



Our integral-field observations uncover complex kine- 
matic structures within the central 3.1 kpc (4.9") of NGC 
6086. In Figure [6j we display 2-dimensional maps of 
kinematic moments from OSIRIS and GMOS data; we 
have computed w ra d, c, h^, and hi by fitting a 4th- 
order Gauss-Hermite polynomial to each non-parametric 
LOSVD. For data with modest signal-to-noise, non- 
parametric LOSVDs must be used with caution, as noise 
may falsely introduce strong non-Gaussian components 
to the fit. In Figure [3 we compare two estimates of i> ra d 
and a from OSIRIS spectra. One estimate is obtained 
by fitting Gaussian profiles to non-parametric LOSVDs 
(Figures 7a and 7b) , and the other is obtained by forcing 
a Gaussian LOSVD to fit the original spectra (Figures 
7c and 7d). In every spatial region, u rac j and a from the 
two fitting options are consistent within errors, and so 
we can trust the non-parametric LOSVDs. GMOS data 
show similar agreement between the two estimates. 

The stellar velocity dispersions measured by OSIRIS 
and GMOS each peak within 250 pc of the galaxy cen- 
ter, but not at the central spatial bin. Central decreases 
in velocity dispersion have been observed in several 
other early-type galaxies with known black holes (e.g., 
van der Marell994lPinknev et al.1l2003t iGebhardt et al.1 
2007t INowak et al.1 12008') . Possible physical explanations 
include an unresolved stellar disk or a localized popula- 
tion of young stars. No dust features are present in pho- 
tometry of NGC 6086, nor is there any evidence of an 
active galactic nucleus. The radial velocities are highly 
disturbed within the central 200 pc, which are only re- 
solve d by OSIRIS: at maxim um, Ai; ra d = 194 ± 52 km 
s _1 . IGebhardt et al.l ()2007f ) found similar patterns in 
u ra d and a in the central 100 pc of NGC 1399, which 
were reproduced by models with a high prevalence of 
tangential orbits. Likcwwise, our best-fitting model of 
NGC 6086 is tangentially biased in the central 200 pc; 
the average value of cr r /a t is 0.55. However, the 2- 
dimensional structure of i> ra d is not consistent with a 
resolved stellar d isk. Axisymme t ric modeling of ellip- 
tical g alaxi es bv IGebhardt et al.l (|2003f) . iShapiro et al.1 
(2006), and lShen fc Gebhardtl (|2010l ) suggests that tan- 
gential bias is common within the black hole sphere of 
influence. 

Figures 6c, 6d, 6g, and 6h illustrate the 2-dimensional 
behavior of the third- and fourth-order Gauss-Hermite 
moments, /13 and /14. Within errors, our measurements 
are largely consistent with /13 = 0, while ft, 4 is signifi- 
cantly negative, corresponding to LOSVDs with "boxy" 
shapes and truncated wings. 

We use major-axis kinematics from ICBH99I to con- 
strain stellar orbit models at radii out to 18.9 kpc (29.3"), 
several times the extent of our IFS data. To incorporate 
these data into our models, we have adopted higher un- 
certainties than the values quoted in ICBH 99; our treat- 
ment attempts to account for additi onal systematic er- 
rors, which are described bv lCBH9"9l but excluded from 
their published measurements for NGC 608 6. In F ig- 
ure we compare kinematic moments from CBH99 to 
the moments derived from OSIRIS and GMOS, selecting 
the spatial bins along the galaxy's major axis. We invert 
the sign of w ra d and /13 for bins on the southern half of 
the galaxy. Values of w ra d, c, and /13 measured from the 
three data sets largely agree, although the OSIRIS data 
yield somewhat smaller values of a. At radii between 
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Fig. 6. — Two-dimensional kinematics in NGC 6086. (a) Radial velocity map from OSIRIS, (b) Velocity dispersion map from OSIRIS, 
(c) Map of ha from OSIRIS, (d) Map of h 4 from OSIRIS, (e) Radial velocity map from CMOS, (f ) Velocity dispersion map from CMOS, 
(g) Map of /13 from CMOS, (h) Map of /14 from CMOS. All maps are derived by fitting 4th-order Gauss-Hermite polynomials to non- 
parametric LOSVDs. The dashed line in each figure represents the major axis of the galaxy, with the arrow pointing north. The solid 
circles in (e) through (h) represent the outermost extent of kinematics from OSIRIS. The median error values are (a) 41 km s" 1 ; (b) 26 
km s" 1 ; (c) 0.042; (d) 0.024; (e) 20 km s' 1 ; (f) 16 km s _1 ; (g) 0.035; and (h) 0.021. 

to the major-axis trends depicted in Figure El and agre e 
with minor-axis kinematics from ILoubser et ail {2008). 
Beyond the central 200 pc, we find no convincing signs 
of kinematically distinct stellar populations dominating 
galaxy spectra at 0.5, 0.9, and 1.6 /xm. 

The most significant discrepancy in the major-axis 
kinematics is between the negative values of /14 derived 
from IF S data and the positive values of /14 measured by 
CBH99 (Figure 8d). Given the uniformity of the positive 
and negative values over many radial bins, we attribute 
this discrepancy to systematic e rrors in at least one set 
of measurements. As noted by CBH99 (and references 
therein), stellar tem plate m ismatch can bias hi] how- 
ever, this study and ICBH99I both perform careful anal- 
ysis with m ultiple stellar temp lates (9 and 28 stars, re- 
spectively) . iNowak et al.l (|2008f) demonstrated that large 
spectral smoothing factors can also bias /14 to negative 
values. Still, a smoothing factor > 100 would be neces- 
sary to produce the full discrepancy between our values 
and those from ICBH99I . Another source of error could 
be adjustments to the equivalent widths of absorption 
features in galaxy and template spectra; we address this 
issue in Appendix [3] Regardless of the cause, includ- 
ing discrepant data in stellar orbit models can influence 
the best-fit solutions. We discuss the effects on measure- 
ments of M. and M+/Lr in Section f4. 3. 21 

4. STELLAR ORBIT MODELS AND BLACK HOLE MASS 

4.1. Stellar Orbits 

We generate stellar orbit models of N GC 6086, using 
the st atic potential method introduced by iSchwarzschildl 
(1979). We use the axisymmetric modeling algo- 
rithm described in detail in Gebhardt et al. (2000b; 
2003) Thomas et al. (2004; 2005), and iSiopis et alj 
( 2009). Here we provide a sum mary of the procedure. 




Fig. 7. — Comparison of Gaussian and non-parametric LOSVD 
extraction methods, (a) Radial velocity map of NGC 6086 from 
OSIRIS, derived by fitting Gaussian profiles to non-parametric 
LOSVDs. (b) Velocity dispersion map from fitting Gaussian pro- 
files to non-parametric LOSVDs. (c) Radial velocity map, derived 
by fitting Gaussian LOSVDs to OSIRIS spectra, (d) Velocity dis- 
persion map from fitting Gaussian LOSVDs to OSIRIS spectra. 
The dashed line in each figure represents the major axis of the 
galaxy, with the arrow pointing north. The median error values 
are (a) 41 km s" 1 ; (b) 26 km s" 1 ; (c) 35 km s" 1 ; and (d) 32 
km s _1 . Radial velocities and dispersions from the two extraction 
methods are consistent within errors. 

0.6" and 4.9", GMOS spectra from the southern half of 
the galaxy yield significantly lower values of a than spec- 
tra from the northern half. The average asymmetry is 
35 km s _1 ; our median error for individual GMOS mea- 
surements is 16 km s _1 . The asymmetry is not seen in 
long-slit data, which more consistently agree with GMOS 
along the north side of the major axis. In spatial bins cor- 
responding to the minor axis, w ra( j and a behave similarly 



Simila r models are prese nt ed in iRichstone fc Tremaind 
(|1984f) . iRix et al.1 (|1997D . iCretton et al.l (|1999t) . and 
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Fig. 8. — Extracted kinematics along the major axis of NGC 6 086. Black circles are OSIRIS measurements, red diamonds are CMOS 
measurements, and blue squares are measurements from ICBH99I . Filled symbols represent points from the positive (north) side of the 
major axis, and open symbols represent points from the negative (south) side. For OSIRIS and CMOS data, the kinematic moments are 
computed from non-parametric LOSVDs. (a) Radial velocity, relative to the central velocity of NGC 6086. The central velocity is defined 
separately for each dataset, such that the spatial averages to either side of the major axis are symmetric about v T& ^ = 0. Velocities from 
the negative (south) side of the major axis have been inverted, (b) Line-of-sight velocity dispersion, (c) Gauss-Hcrmite /13, with values 
from the negative (south) side of the major axis inverted, (d) Gauss-Hermite /14. 



iValluri. Merritt fc Emselleml (|2004D . 

We assume that the central region of NGC 6086 con- 
sists of three mass components - stars, a central black 
hole, and an extended dark matter halo - described by 
the radial density profile 



closed halo mass for this profile, 



M halo (< r) 



G 



1 



(4) 



p(r) 



Mi 

Lr 



v{r) + M.5(r) + phaioW . 



(2) 



The stellar distribution is assumed to follow the observed 
(deprojected) luminosity density v{r) (see Figure [IJ 
with a constant stellar mass-to- light ratio M^/L^. For 
the dark matter halo, we com pare two density profiles: 
the commonly used NFW form (jNavarro. Frenk fe; White! 
fl996t) . and a logarithmic (LOG) profileQ: 



PhaloM 



3ri 



4irG ( r 2 + r 2y 



(3) 



The free parameters in the LOG profile are the asymp- 
totic circular speed v c and the core radius r c , within 
which the density is approximately constant. The en- 



The gravitational potential is logarithmic: Q 



ln(r 2 



is predominantly set by v c . The difference between the 
NFW and LOG profiles is greatest at small radii, where 
the NFW profile yields higher densities, phaio cx r -1 . 
However, each profile is greatly exceeded by the stellar 
mass density in the inner regions of NGC 6086. As de- 
scribed below, we have compared LOG and NFW profiles 
in a subset of models, and find no significant differences 
in the best- fit values of M # . 

For a given set of input parameters M # , M+/Lf>, and 
Phaio, we compute a continuous, static gravitational po- 
tential from Equation Azimuthal symmetry about 
the z-axis (corresponding to the projected minor axis) is 
imposed, as well as symmetry about the equatorial plane 
(z = 0). We then generate stellar orbits by propagating 
test particles through the potential. Orbits are tracked 
in a finely spaced polar grid, {r,6), where 9 is the po- 
lar angle from the z-axis. Our models of NGC 6086 use 
96 radial and 20 polar bins per quadrant. Each orbit is 
sampled at a random set of azimuthal angles, 4>. 

The initial phase space coordinates of test particles 
are chosen to sample thoroughly three integrals of mo- 
tion: energy E, angular momentum component L z , and 
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the third, non-classical integral, I3 . Computational noise 
and finite propagation steps introduce noise into test par- 
ticle trajectories; this is mitigated by allowing each par- 
ticle to complete 200 circuits of the potential and then 
determining its average orbit. Orbits that escape the 
potential are not included in subsequent fitting. For a 
given potential, our model of NGC 6086 includes ap- 
proximately 16, 000 to 19, 000 bound orbits. Identical 
counterparts with the opposite sign of L z raise the to- 
tal to 32,000 - 38,000 orbits. Each orbit in the model 
is assigned a scalar weight; initially, all bound orbits are 
given equal weights. 

The set of best-fit orbital weights is determined by 
comparing projected LOSVDs from the orbits to the ob- 
served LOSVDs for the galaxy. Our models use non- 
parametric LOSVDs, defined in 15 velocity bins between 
-1000 and 1000 km s" 1 . Each observed LOSVD spa- 
tially maps to a linear combination of bins within the 
model, according to the spatial boundaries of the corre- 
sponding spectrum, and to the instrument-specific PSF. 
A corresponding model LOSVD is computed from the 
projected velocity distributions of individual orbits in 
each spatial bin, the appropriate combination of spatial 
bins, and the orbital weights. Only the orbital weights 
are varied to determine the best-fit solution. 

The best-fit solution is determined by the me t hod o f 
maximum entropy, as in iRichstone fe Tremaind ((1988). 
This method maximizes the function / = S — ax 2 , where 

9 [^i.data ( v j) ~ C>i, model i v j)] /r \ 

* =2^2^— _ 2 ,.. \ — (5) 



» j 



°f(«j) 



and 



Here, C 



S = - w k In 



(6) 



z,data 



and £ 



i, model 



are LOSVDs in each of the i = 
1, ... Nf, spatial bins, o~ 2 (vj) is the squared uncertainty in 
A, data at velocity bin Vj, Wk is the orbital weight for the 
fcth orbit, and Vk is the phase volume of the fcth orbit. 
The parameter a is initially small so as to distribute 
orbital weights broadly over phase space, and is increased 
over successive iterations so that the final optimization 
steps exclusively minimize x 2 . A further constraint for 
all solutions is that the summed spatial distribution of 
all weighted orbits must match the observed luminosity 
density profile. 

4.2. Black Hole Mass and Mass-to-Light Ratio 

Our combination of IFS and long-slit data within 30" 
is not sufficient to measure directly the dark matter halo 
profile or enclosed mass. We therefore consider different 
halo masses and profiles in our analysis. For a given dark 
matter halo profile, we generate a set of 32,000-38,000 
stellar orbits for an input M, and M*/L#, and obtain 
X 2 for the best-fit orbital weights using the method de- 
scribed in the previous subsection. This process is re- 
peated over a finely sampled grid in M. and M*/Lr. We 
have completed all trials for NGC 6086 on supercomput- 
ers at the Texas Advanced Computing Center (TACC), 
totaling ~ 10, 000 CPU hours. Our results are summa- 
rized in Tabled 



In Figure \§\ we illustrate how x 2 i n the orbit mod- 
els varies with M, and M^/L^. Three dark matter halo 
masses are shown: no dark matter (left panel), an inter- 
mediate LOG halo with v c = 300 km s _1 (middle), and 
a maximal LOG halo with v c = 500 km s _1 (right). The 
latter is chosen to approximate the measured line-of-sight 
velocity dispersion of 302 km s _1 for NGC 6086's host 
cluster (|Zabludoff et al.l fl993f ). which corresponds to a 
full 3-dimensional velocity dispersion of 523 km s _1 . We 
set the core radius in Equation © to be r c = 8.0 kpc, re- 
flectin g the value of 8.2 kpc determined by Thomas ct al. 
(2007) for NGC 4889, the Coma BCG. Our dynamical 
models are constrained within a radius of 18.9 kpc, cor- 
re sponding to outermost radius of 29.3" for long-slit data 
in CBH99. We also have run models with a single NFW 
dark matter profile, constructed to contain the same en- 
closed mass within 18.9 kpc as our most massive LOG 
halo. In cosmological A-body simulations, the NFW 
scaling parameters c and r s are correlated according to 
the relationship 



r 3 = 



3 x 10 13 M 



o 



200f p crit c 3 



10o35 (1-05-logioc) 



(7) 



((Navarro. Frenk fe Whitil996l: iRireTalinMFl where 
pah — 3Hq/8ttG. Combining this relationship with our 
enclosed mass constraint, we obtain c = 9.6 and r s = 94.0 
kpc. 

Dark matter is ubiquitous in galaxies, thus motivat- 
ing its inclusion in stellar orbit models. Furthermore, 
models with dark matter produce better fits to our full 
set of kinematics: when the dark matter component is 
removed, Xmin increases by ~ 100 (Table [2]). In Ap- 
pendix |B] we compare our full sets of observed LOSVDs 
to the best-fitting models with and without dark matter. 
The largest discrepancies between the two models occur 
at ra dii beyon d 15", where we only have a few data points 
from lCBH99l . Without thorough radial coverage or mul- 
tiple long-slit position angles, we cannot fully untangle 
degeneracies between M+/Lr, v c , and r c (or M*/Lr, c, 
and r s in the case of an NFW profile) . 

In Figure [Til] we display the dark matter fraction as 
a function of radius, for each of the halos described 
above. In each case, we use the best-fit values of M, 
and M+/Lr, described below, to compute the total en- 
closed mass. Using our surface brightness profile from 
HST /KPNO we compute an effective radius, R c g, of 
31.7" (20.4 kpc), defined as the semi-major axis of the 
elliptical isophote containing half of the total luminosity. 
We assume a tot al luminosity of 1.82 x 10 11 L© R, from 
M v = -23.11 in lLauer et all ((20071 ) and V — R = 0.64. 
Within Reg, dark matter composes 44% to 74% of the 
total mass (for the intermediate-mass LOG halo and the 
NFW halo, respectively). This range agrees with dynam- 
ical m odels of other c D gal axies with LOG and NFW 
halos: iThomas et all ((20071 ) found - 50 - 75% dark 
matter within R e g for NGC 4889 and NGC 4874, and 
IGebhardt fc Thomasl J2009) found - 40% dark matter 
within R eS for M87. Within the 4.9" outer radius of 
GMOS data, the maximum dark matter fraction in our 
models is 20%. 



11 Note the erratum in Equation (B3) of IHix et all (1997); the 
correct equation is logioc = 1.05 — 0.15 logio (M200/3 X 1O 13 M0). 
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Fig. 9. — Surface plots of x 2 vs. M*/Lr and M., for models fitting both OSIRIS and CMOS data. Left: Models with no dark matter 
halo. Long vertical pixels represent regions with coarser sampling in M m . Middle: v c = 300 km s _1 and r c = 8.0 kpc, for Mh a i = 3.3 X 10 11 
Mq within 18.9 kpc. Right: v c = 500 km s _1 and r c = 8.0 kpc, for M halo = 9.3 X 10 11 Mq within 18.9 kpc. For each dark matter halo, 
additional models were run outside the range of M^/L^j and M m depicted here. These models all yield higher values of \ 2 . 



TABLE 2 

Results from stellar, orbit models. 





Data 




(km s- 1 ) 


(kpc) 




(1) 




(2) 


(3) 


OSIRIS - 


h GMOS - 


h CBH99 


0.01 


20.0 


OSIRIS - 


b GMOS - 


h CBH99 


300 


8.0 


OSIRIS - 


b GMOS - 


h CBH99 


500 


8.0 


OSIRIS - 


b GMOS - 


h CBH99 







c r s 
(kpc) 

(4) (5) 



iWhalo 
(10 11 Mq) 

(6) 



(10 9 Ms) 
(7) 



M*/L fl 
(M /L ) 

(8) 



Xmin 



(9) (10) 



OSIRIS + GMOS 
OSIRIS + GMOS 
OSIRIS + CBH99 
GMOS + .CBH99 



0.01 
500 
500 
500 



20.0 
8.0 
8.0 
8.0 



9.6 94.0 



« 
3.34 
9.28 
9.05 
« 
9.28 
9.28 
9.28 



0.6Z 



2.6 



+1.0 
l.O 



3.6Z 



3.6^ 



3.2 



+1.5 
1.3 



3.5: 



1.9- 1 



0.6Z 



3.6 



hi. 3 

-1.0 



+ 1-T 



3.6 



+ 1.6 



3.2 



+ 1.8 
1.3 



3.5 



+ 1.8 



+ 1.5 



1.9 



6.68 
5.1 
4.6 
4.3 
4.4 
4.lj 
4.6^ 
4.5 4 



+0.15 
0.17 



..7 



+0.2 



-0.3 
-0.3 



-0.3 
-0.4 



-0.3 
-0.4 



+0.4 
0.3 



s +0.7 
J -0.7 



5.4 
14.6 
^4.3 
^4.4 
(4.1 
^4.6'_ 
f4.5^ 



+0.3 
0.8 



+0.3 \ 
0.7 J 



+0.3 \ 
0.7 i 



0.6 ) 



+0.A\ 
0.7 i 



+0.7> 



1114.4 795 

1019.4 795 

1010.2 795 

1009.4 795 

622.9 345 

623.8 345 

412.7 570 

725.5 675 



Notes: Column 1: Data sets included in trial. Column 2: circular velocity of LOG dar matter halo (Eq. |3}- Column 3: 
core radius of LOG dark matter halo (Eq. [3]). Column 4: concentration parameter for NFW dark matter halo. Column 5: 
scale radius for NFW dark matter halo. Column 6: enclosed halo mass, defined at the outermost long-slit data point. The 
corresponding radius is 18.9 kpc. The "no dark matter" case has v c = 0.01 km s _1 , r c — 20.0 kpc, and Mhalo ~ 200 Mq. Column 
7: best-fit black hole mass. Quoted errors correspond to 68% confidence intervals. Values in parentheses include all systematic 
errors. Column 8: best-fit i?-band stellar mass-to-light ratio. Quoted errors correspond to 68% confidence intervals. Values in 
parentheses include all systematic errors. Column 9: minimum x 2 value for all models. Column 10: degrees of freedom in model 
fits to LOSVDs. Computed values include a smoothing factor of 1 degree of freedom per 2 velocity bins for non-parametric 
LOSVDs from OSIRIS and GMOS. 
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Fig. 10. — Dark matter fraction of total enclosed mass, for the 
model dark matter halos presented in Table [2] For each halo, the 
total enclosed mass includes the best-fit black hole mass and best- 
fit stellar mass-to-light ratio, derived using OSIRIS, CMOS, and 
long-slit data. The dashed vertical line marks the effective radius 
derived from our photometric data. Filled symbols represent the 
radi al positi ons of spectra from OSIRIS (circles), CMOS (stars), 
and CBH99 (squares); the vertical positions of the symbols have 
no physical meaning. 

In Figure [TT] we illustrate the variation of M, and 
M^/Lr with enclosed halo mass. We compute best- 
fit values for M, and M^/L^ by integrating the 2- 
dimensional likelihood function from each % 2 surface; we 
describe this method and our determination of errors in 
Section 14.31 Figures |H] and [11] show that the best- fit val- 
ues of M, and M^/Lr are substantially influenced by the 
presence of dark matter in the stellar orbit models. This 
occurs because our innermost kinematics sample the nu- 
cleus of the galaxy, where orbits are dominated by the 
enclosed mass of stars and the central black hole, whereas 
enclosed stellar and dark halo masses are both important 
at larger radii. The significant presence of dark matter 
at large radii drives the best-fit models to lower values 
of M*/L#. In turn, the decreased stellar mass requires 
a higher black hole mass to reproduce the kinematics in 
the nucleus. This tr e nd wa s initially demonstrated by 
IGebhardt fc Thomas! ((2009), for M87. Ne gative covari- 
ance between M. and M^/L^ is also visible in % 2 con- 
tours for individual dark matter halo models (Figure [9]). 
In Section [SJ we compare our best- fit values of M, from 
different dark matter halo models to the predictions from 
the M, — a and M, — L relationships. 

Using an NFW profile yields a 7% decrease in the best- 
fit value of M*/Li{, relative to the LOG profile with the 
same enclosed mass at 18.9 kpc. This is because the 
majority of our kinematic measurements occur at r < 5 
kpc, where the more centrally-concentrated NFW profile 
yields higher dark matter densities. In the central 100 pc, 
the stellar core of NGC 6086 varies nearly as r _1 in lumi- 
nosity density, mimicking the slope of the NFW profile. 
Near the black hole, the lower stellar mass density bal- 
ances the higher density in dark matter, and so the best- 
fit values of M, are identical for the NFW and LOG pro- 
files. We find X 2 mi n,WG ~ Xmin,NFW = °- 8 ; indicating no 
significant d iffer ence in the goodness of fit. iThomas et all 
(|2005l[2007l) and lGebhardt fc Thomasl (|2009D have found 
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Fig. 11. — Top: Best-fit black hole mass for different assumed 
dark matter halos. Bottom: Best-fit stellar mass-to-light ratio. 
Black circles represent LOG halos, and grey diamonds represent 
the NFW halo. Each measurement co mes from fitting the full set 
of LOSVDs from OSIRIS, CMOS, and lCBH99l The enclosed halo 
mass is define d at 18.9 kpc, the radius of the outermost data point 
from ICBH99I . Error bars indicate 68% confidence, and include 
systematic errors. 

similar difficulties in distinguishing between NFW and 
LOG profiles. 

A second way to address the influences of dark mat- 
ter on M+/Lr and M, is to fit the orbit models only at 
radii where dark matter composes a small fraction of the 
enclosed mass. For NGC 6086, we have run two trials 
in which we only fit LOSVDs from OSIRIS and CMOS: 
one trial with the maximum-mass LOG dark matter halo 
described above, and one trial with no dark matter. In 
both of these trials, the best-fit values of M. and M^/Ln 
agree with our results from fitting IFS and long-slit data 
with the maximum-mass LOG halo (see Table [2]). This 
agreement provides strong evidence that M. ~ 3 x 10 9 
M Q , regardless of our insensitivity to the exact structure 
of the dark matter halo in NGC 6086. In contrast, forcing 
models without dark matter to fit long-slit data biases 
the best-fit black hole mass to a substantially lower value, 
~6x 10 8 Mq, and increases the best-fit stellar mass-to- 
light ratio to 6.7 Mq/L q . This mass-to-light ratio is 
highly inconsistent with stellar population estimates, as 
we discuss in Section [5] Even though excluding long- 
slit data results in consistency between models with and 
without dark matter, these models are not as thoroughly 
constrained, and we obtain slightly larger confidence in- 
tervals in M, and rVL/L^ for each trial 

Our data from OSIRIS, CMOS, and |Clffl99| play com- 
plementary roles in constraining the gravitational poten- 
tial of NGC 6086. In Figure [12] we compare model re- 
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suits for one dark matter halo (LO G; v c = 500 km s 1 ), 
using data only fr om GM OS and ICBH99I . versus only 
from OSIRIS and ICBH99I . The GMOS data are suffi- 
cient to detect a black hole, in part because of excellent 
seeing. Yet the strong diagonal contours in the left panel 
of Figure indicate that the black hole mass derived 
from GMOS is degenerate with the enclosed stellar mass. 
LOSVDs from OSIRIS have large statistical errors and by 
themselves cannot place strong constraints on the black 
hole mass. However, the OSIRIS data help separate the 
respective influences of the stars and the black hole. Us- 
ing GMOS and OSIRIS data together reduces covariance 
between M, and M*/L/j and lowers the statistic al uncer- 
tainties of both quantities. Long-slit data from [CBH99 
confirm the presence of dark matter and tighten con- 
straints on M*/Lft and M m for individual dark matter 
halo models. 

It is not clear how to interpret the significant increase 
in the best-fit value of M„ from x 10 9 M Q with 

OSIRIS and long-slit data only, to 7±| x 10 9 M with 
GMOS a nd long-slit data only. Using central a = 329 km 
s- 1 from lLoubser et all (|2008l) . and M. = 3.6 x 10 9 M 
from our combined-data trial, we compute r; n f = 0.22". 
In this case, GMOS marginally resolves the sphere of in- 
fluence, with a seeing FWHM ~ 2r in f. At small radii, 
LOSVDs from GMOS yield slightly higher velocity dis- 
persions than overlapping LOSVDs from OSIRIS (Figure 
8b) , which could contribute to the increase in M, . With 
no obvious way to assess independently the accuracy of 
data from OSIRIS versus GMOS, we favor including both 
sets of LOSVDs. The corresponding black hole mass is 
3.6 jli' i x 10 9 Mq, which lies between the two partial- 
data values and has the narrowest confidence interval. 
The confidence interval in M*/L# is also minimized by 
including all data, though the best-fit value of 4.6l ; 7 
does not change significantly upon exclusion of OSIRIS 
or GMOS data. 

None of the \ 2 surfaces in Figures |H] and [T^] are 
completely smooth. In particular, the models without 
any dark matter show large variations in x 2 over small 
changes in M, and M^/L^ (Figure [9]). We suspect that 
these variations arise from numerical noise in propagat- 
ing test particles through different potentials: each small 
change in the potential may send a given test particle 
through a different set of spatial regions. The cumula- 
tive effect is that each model creates a different set of 
test-particle LOSVDs, and x 2 can change abruptly in 
spite of the freedom to adjust orbital weights. When 
the models include a constant dark matter component, 
the relative changes in the potential are smaller, and the 
noise in x 2 is less pronounced. Still, the x 2 surface for 
each dark matter halo exhibits a noise floor at the level 
of Ax 2 ~ 1. In Section l4~3Jl we describe how this noise 
influences our measurements of confidence intervals. 

4.3. Determining Errors 

The figure of merit for evaluating confidence intervals 
in M. and M*/Lr is Ax 2 = X 2 - Xmi n > where Xmi n is 
the lowest output value among all models. For NGC 
6086, we determined confidence intervals by integrating 
the relative likelihood function, P a e"i Ax . Although 
Ax 2 is a bette r statistical indicator than y 2 per degree of 
freedom (e.g.. Ivan der Marel et al.lll~99 8; Gcbhardt et al.l 



2003), the latter is useful for crudely indicating the level 
of agreement between the data and the model with the 
best fit. For each model, the number of degrees of free- 
dom, Adof , depends on the number of observed LOSVDs 
and the number of velocity bins evaluated per LOSVD. 
Because we used a spectral smoothing factor in deter- 
mining LOSVDs for OSIRIS and GMOS, the velocity 
bins are not entirely independent. We estimate that each 
LOSVD from OSIRIS or GMOS has 1 degree of freedom 
per 2 velocity bins, whereas long-slit data from[CBH99 
has 1 degree of freedom per velocity bin. For all exper- 
iments with NGC 6086, we find Xmin/^dof between 0.7 
and 1.8, indicating reasonable agreement. 

4.3.1. Confidence Intervals 

We determine confidence intervals by using Ax 2 as an 
empirical measure of relative likelihood between models 
with different M, and M t /Ljj, and by numerically inte- 
grating this likelihood with respect to M, and M*/L^. 
In contrast to the majority of previous studies (e.g., Geb- 
hardt et al. 2000b; 2003; 2007; Nowak et al. 2007; 2008; 
Gebhardt & Thomas 2009; Giiltekin et al. 2009b; Siopis 
et al. 2009; Shen & Gebhardt 2010; cf. van der Marel 
et al. 1998), we do not use fixed values of Ax 2 to define 
confidence intervals. The fixed Ax 2 method is appropri- 
ate only if the orbit models cleanly sample a well-defined 
likelihood function of M. and M^/L^; other studies typ- 
ically assume a 2-dimensional Gaussian likelihood func- 
tion. Our models of NGC 6086 produce noisy x 2 con- 
tours (Figure [9]), in which case the fixed Ax 2 method 
is sensitive to noise in individual models. This effect is 
especially pronounced for models without a dark matter 
halo. 

We define likelihood P such that two models with x 2 
and x| have relative likelihood 



Pi 

P-2 



(8) 



This form of P is valid, provided that x 2 is measured 
from indepen dent, Gaussian-distributed data points 
(|Cowanl [[998); in our case, these are the observed 
LOSVDs. To evaluate likelihood with respect to a single 
variable (i.e. x = M.), we marginalize the 2-dimensional 
surface with respect to the other variable (y = M*/L,r), 
such that: 



Krnax 

P(x) oc e^* 2 ^ Sy , (9) 

J/min 

where Sy is the interval between sampled values of y. 
Confidence intervals in x are determined by evaluating 
the cumulative distribution: 



C(x) 



f^ iTi P(x')dx> 



(10) 



In practice, we define [x mln , x max ] and [y min , y max ] by ex- 
panding our range of models until the marginalized like- 
lihood functions P(M.) and P(M*/L^) are nearly zero 
at the minimum and maximum modeled values of M. 
and M^/L/j;. The physical limit M. = is included in all 
trials. For confidence level fc, we define confidence lim- 
its at C = | (1 ± k). For example, the 68% confidence 
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Fig. 12, — Surfa ce plots of \ 2 vs - M*/Lfl and M., using integral-field data from different instruments. Left: GMOS data plus long-slit 
data from CBH99| . Right: OSIRIS data plus long-slit data from CBH99. All models include a LOG dark matter halo with v c = 500 km 
s~ 1 and r c = 8.0 kpc, yielding Afi la i Q = 9.3 X 10 11 Mq within 18.9 kpc. In each case, additional models were run outside the range of 
M t /Lij and Mm depicted here. These models all yield higher values of \ 2 ■ 



interval comprises all x for which 0.16 < C < 0.84. For 
each trial, we define the best-fit values of M. and M*/L/j 
as the median values from P(M,) and P(M*/Lr), cor- 
responding to C = |. In Figure HH we show P(M,) 
and P(M- k /Ln) for each dark halo setting, along with 
cumulative distributions, median values, and confidence 
intervals. To estimate precise confidence limits, we lin- 
early interpolate C between discretely sampled values of 
M m and M^/Lr. 

Our empirical treatment yields wider 68% confidence 
intervals than those derived from fixed A\ 2 and a Gaus- 
sian likelihood function. Our intervals for confidence lev- 
els > 90% typically fall near those derived from fixed 
A% 2 . By construction, our confidence intervals do not 
include M, = 0. For the maximum-mass LOG dark 
matter halo in NGC 6086, the marginalized likelihood 
corresponding to M. = 0, P(M. = 0), is 0.06% of the 
maximum marginalized likelihood value. For a Gaussian 
likelihood function, this likelihood ratio would indicate 
> 99.98% confidence for a black hole detection. For mod- 
els without dark matter, P{M m = 0) is 23.2% of the max- 
imum value, and the detection falls to 91% confidence. 

An alternative way to determine confidence intervals 
from noisy \ 2 is the method of Ivan der Marel et al.l 
(1998), in which random noise is added to the LOSVDs 
output by each orbit model, the \ 2 surface is re- 
computed, and confidence limits arc determined using 
fixed Ax 2 . This process is repeated in a Monte Carlo 
fashion, and the extrema of the confidence limits from all 
trials are adopted. This treatment assumes that numeric 
noise in orbit models produces fluctuations about an in- 
trinsically Gaussian likelihood distribution. The global 
likelihood function of our models, however, is visibly non- 
Gaussian (see Figure H3l top, in particular). 

4.3.2. Systematic Errors 

The confidence intervals measured from x 2 surfaces ac- 
count for statistical errors in the observed LOSVDs and 
random noise within the stellar orbit models. Systematic 
errors must be addressed separately. We have directly 
tested several systematic effects. 



Our largest systematic error arises from discrepancies 
between LOSVD shapes derived from IFS versus long-slit 
data, as indicated by the parameter /14 (Figure 8d). To 
test how this discrepancy biases M, and M+/L r, we fit 
Gaussian profiles to the orig inal LOSVDs from CBH99, 
constructing an alternative set of LOSVDs with /13 = 
and hi = 0. Using these new LOSVDs in combination 
with the OSIRIS and GMOS data, we repeated the stel- 
lar orbit models for the maximum dark matter halo, and 
found Mm = 4.5±^_f x 10 9 Mq, and M±/L R = 4.2±°j 
M e /L0. This is a 22% increase in M, and a 9% de- 
crease in M^/Lfl, relat ive to th e corresponding trial with 
original LOSVDs from lCBH99l . The direction of the bias 
indicates that more enclosed mass is required to produce 
LOSVDs with extended wings (h± > 0) at - 3 - 19 
kpc. Since the additional mass is obtained by increasing 
M*/Lfl, the innermost LOSVDs drive the best fit toward 
a smaller value of M, . 

Two additional sources of systematic error are the 
uncertainties of the optimal stellar template for fitting 
LOSVDs, and the average PSF for OSIRIS. We find that 
both effects yield small errors. Trials with 3 PSFs yield 
a standard deviation of 8% in M, and 1% in M*/Lfl. 
Trials with 2 stellar templates differ by 5% in M, (corre- 
sponding to 3.6% standard deviation) and 2% in M^/L/j 
(1.6% standard deviation). We describe these tests in 
detail in Appendix [Al 

We use the following prescription to compute the total 
error in M, for a given dark halo: 



c+, tot 



'PSF 



^temp + 2< 5/i 4 



and 



-, tot 



= « 



°PSF 



"'temp) 2 



(11) 

(12) 



Here, <j+ t tot and o- t tot are the upper and lower portions 
of the 68% confidence interval, including all errors; er + x 2 
and ct_ i x 2 are the contributions from integrating the em- 
pirical likelihood function as in Section |4. 3. 1[ and repre- 
sent statistical errors; itpsf is the standard deviation in 
best-fit Mm from trails with different OSIRIS PSFs; a tcmp 
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Fig. 13. — Likelihood functions P(M,) and P(M*/L.r), after marginalizing over one variable (M*/Lfl and M., respectively) and re- 
normalizing. The solid red line in each plot is the cumulative distribution. Vertical lines in each figure represent the median value (dotted 
red line) and intervals for 68%, 90%, and 99% confidence (dashed lines). Top: Models with no dark matter halo. Middle: v c = 300 km 
s" 1 and r c = 8.0 kpc, for M halo = 3.3 X 10 11 M Q within 18.9 kpc. Bottom: v c = 500 km s^ 1 and r c = 8.0 kpc, for A/ halo = 9.3 X 10 11 
M Q within 18.9 kpc. 

is the standard deviation from trials with different tem- 
plate stars; and 5h 4 is the difference in best-fit M , from 
using re-fit versus original LOSVDs from CBH991 As ex- 
cluding I13 and hi from these LOSVDs introduces a bias 
toward higher M., we assign this effect solely to <r+,tot, 
with a magnitude of 2 x "^^ 4 - Our equations defining 

er +! tot and cr_ ! to t for M^/Lr are similar to (TTTI) and (|T2|) . 
except we apply dh t entirely to a^.tot- To apply our re- 
sults to different dark matter halo settings, we define the 
systematic errors as percentages, such that Sh 4 , opsf and 
Ctemp scale with M, and M^/L^. We list the 68% con- 
fidence intervals, including all errors, inside parentheses 
in Table [2] The confidence intervals outside parentheses 



in Table [2] only include cr^. x2 . 

By adding in quadrature, we have assumed zero corre- 
lation between different sources of systematic error; this 
is the most conservative appro ximation . <5/ t4 and Ctcmp 
are likely correlated (see, e.g., ICBH99I ). but the contri- 
butions from (Ttemp are small, and we have not run ex- 
tensive tests to measure covariance between stellar tem- 
plates and overall trends in h^. There is no obvious rea- 
son to suspect covariance between other terms. 



4.3.3. Other Potential Sources of Error 



The shapes and depths of CO bandheads depend 
on luminosity class a s well as spectral type (e.g., 
iSilge fc Gebhardtl 12003). Although we have extensively 
examined stars of multiple spectral types, our library 
of template stars from OSIRIS lacks an M-dwarf tem- 
plate; it is uncertain whether dwarf stars contribute sig- 
nificantly to spectra of NGC 6086. Our sample of tem- 
plates is also limited by the range of elemental abun- 
dances found i n brig ht stars in the solar neighborhood. 
lLoubser et al.l ([2009D determine stellar metallicity at the 
center of NGC 6086 to be - 2 times solar ([Z/H] = 
0.28 ± 0.07), with an a-enhancement ratio, [E/Fe], of 
0.39 ± 0.04. Other BCGs are similarly m etal- and a-rich 
(IBrough et al]l2007tlLlivibser et al.l2009l) . If LOSVDs de- 
rived from Call, CO, and OH absorption features are sen- 
sitive to template star metallicities and a ratios, then our 
systematic error could be higher than estimated above. 

Another issue is the possibility of significant spatial 
variations in the stellar mass-to-light ratio, contradicting 
the uniformity imposed upon stellar obit models. Radial 



been measured in individual BCGs 1 


BrouEh et al.l 120071 


iCoccato. Gerhard & Arnaboldi I201C 


); further modeling 



is necessary to quantify corresponding gradients in M*/L. 

An untested source of systematic error for NGC 
6086 is the assumed shape and inclination of the 
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galaxy. Our orbit models use an oblate, axisym- 
metric potential and assume edge-on inclination. 
Previous studies have indicated that uncertain in- 
clination can bias M, by 30 — 50% in elli ptical 
galaxies (iVerolme et all 120021 : iGebhardt etall (20031 : 
IShapiro et al.l 120061 ). Moreover, isophotal evidence and 
simulations of galaxy merg ers suggest that many BCGs 
are prolate or triaxial dPorter. Schneider fc Hoessell 
19911: iRvden, Lauer fc Postman! 119931: 

Boylan-Kolchin. Ma fc Quataertl 120061 ). Although 
the orbit superposition method was originally devel- 
oped for a triaxial potential (Sch warzschildl [1979), 
triaxial orbit models with the spatial and velocity 
resolutions necessa ry to measure M, are a very re- 
cent development (|van den Bosch et al.l [20081. In an 
early comparison of triaxia l and axisymmetric models, 
Ivan den Bosch fc de Zeeuwl (|2010f ) found that M, was 
unchanged in M32, whereas triaxial models of NGC 
3379 increased M, by a factor of 2. 

5. CONCLUSIONS AND DISCUSSION 

We have reported the first stellar dynamical measure- 
ment of the central black hole mass in a BCG beyond the 
Virgo Cluster. Our results are based on 2-dimensional 
stellar kinematics in the central region of NGC 6086 (the 
BCG in Abell 2162): the inner 0.9" (580 pc) from the IFS 
OSIRIS with LGS-AO at Keck, the inner 4.9" from the 
IFS GMOS- N at Ge mini North, and long-slit data out 
to 30" from ICBH991 The individual datasets play com- 
plementary roles in constraining the gravitational poten- 
tial of NGC 6086. Used together, CMOS and OSIRIS 
data reduce degeneracy between the black hole mass and 
enclosed stellar mass near the center, decreasing the un- 
certainties of both quantities. The long-slit data confirm 
the presence of dark matter and constrain the total en- 
closed mass, but they are insensitive to the precise form 
of the dark matter halo profile. 

We have used axisymmetric stellar orbit models in- 
cluding a dark matter halo to determine M m and the R- 
band stellar mass-to-light ratio. We have tested several 
dark matter halo profiles with our full set of kinematic 
measurements; in each case, the best-fit black hole is 
at least 4 times as massive as the best fit without dark 
matter. Including dark matter in the models decreases 
the best-fit value of M+/L R by 20 - 40%. For the most 
massive halo allowed within the gravitational potential 
of the host cluster, we find M. = 3.6t^ x 10 9 M and 
M*/L_r = 4.6±Q y Mq/Lq. We obtain similar values of 
Mm and M*/L# when we exclude long-slit data, for mod- 
els with and without dark matter. 

In Figure HH we add our meas urement of M, in NGC 
6086 to the M.-a relationship of lGiiltekin et al.l (|2009t ). 
Our plotted measurement corresponds to the most mas- 
sive dark matter halo used in our models. We derive 
an effective velocity dispersion of 318 ± 2 km s _1 within 
1 effective r adius, by weighting the kinematic measure- 
ments from ICBH99I with respect to our measured sur- 
face brightne ss profile. The M, — a relation and intrinsic 
scatter from IGiiltekin et al.l (|2009l ) yield a er-predicted 
black hole mass of 0.9^J'g x 10 9 M©. Our measurement 

of Mm = 3.6lti" i x 10 9 M Q with a maximum-mass dark 
matter halo is marginally consistent with this prediction, 
while our measurement of M. = 2.6 jli o x 1O 9 M for an 
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Fig. 14. — The M, — a relationship. The filled circle is 
our measurement of M, = 3. 6±{^ x 10 9 M in NGC 6086, 
using the maximum-mass LOG dark matter halo (n c = 500 
km s — 1 ) and our full set of data from OSIRIS, GMOS, and 
IC-BH991 Without dark matter, we measure Mm = 0.6±°;t x 

10 9 Mq. The remaining data points are the measurements 
compiled in Giiltckin ct al. (200 3), plus updated measurements 
for M 87 and M60 l IGebhardt fc' Thomas 2009; Shcn & Gebhardti 
[201Cf) . The thick dashed line is the fit log(M.) = 8.12 + 
4.24 log (ct/200 kms -1 ) , and the dotted lines enclose a scatter of 
0.44 dex, from lGiiltekin et al.1 | |2009T) 

intermediate-mass halo is closer to the predicted value. 
Our measurement of M. = 0.6 ± 0.4 x 10 9 M without 
dark matter also agrees with the predicted value. How- 
ever, the unrealistic omission of dark matter at radii cov- 
ered by long-slit data biases M, toward low values. The 
V-band luminosity of NGC 6068 is 1.4 x lO 11 L Q v , us- 
ing My = -23.11 from lLauer et~ap (12001. The M. - L 
relation and intrinsic scatter from lGiiltekin etall (p009w 
yield a prediction of M, = 1.3^q' 9 x 10 9 M , which is 
consistent with our measurements with and without dark 
matter. 

Of the existing sample of ~ 30 galaxies with 
Aim measured from stellar dynamics, this study of 
NGC 6086 is only the third to consider dark mat- 
ter. Like NGC 6086, measurements of M, and M*/L 
in M87 depe nd strongly on the inclusio n of a dark 
matter halo (|Gebhardt fc Thomasl 120091 ). However, 
iShen fc Gebhardtl (|2010f) found that introducing dark 
matter to models of M60 produced minimal changes in 
Mm and M*/L. Several factors contribute to the greater 
importance of dark matter in models of NGC 6086 and 
M87. First, the quality and spatial resolution of the 
kinematics are insufficient to distinguish thoroughly the 
black hole from the central stellar mass profile. Neither 
NGC 6086 nor M87 have good spectra from HST, and 
the potential advantage of AO IFS data for NGC 6086 
was compromised by low S/N and coarse spatial binning. 
In principle, very high-quality data could reveal unam- 
biguously a compact mass at the center of the gravita- 
tional potential, permitting robust measurements of M. 
in spite of biased M*/L# values. Second, shallow stel- 
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lar mass profiles in NGC 6086 and M87 limit the range 
of radii where stars dominate the gravitational potential, 
especially in the presence of a relatively massive and con- 
centrated dark matter halo. In contrast, models of M60 
can exclude dark matter and still accurately constrain 
M*/L by removing only the few outermost data points. 
Third, shallow stellar light profiles affect kinematic mea- 
surements near the centers of NGC 6086 and M87: there 
is more contamination along the line of sight from stars 
at larger radii, where dark matter is influential. The lat- 
ter two factors suggest that existing stellar dynamical 
measurements of M, are most likely biased in galaxies 
with large stellar cores. New and revised measurements 
using models with dark matter could raise or steepen the 
upper end of the black hole scaling relations, as stella r 
core size increases with galaxy mass (|Lauer et alJ l2007). 

Stellar population modeling can yield an independent 
measurement of the stellar mass-to-light ratio. To com- 
pare existing stellar population studies to our dynami- 
cal results, we translate M*/Ly and M*/Lj to R band 
by using the {g — r) and (r — i) colors of NGC 6086 
from the Sloan Digit al Sky Survey (SPS S '). App lying the 
filter translations of IBlanton fc Roweisl (|2007f) . we find 
V — R = 0.64 and R - I = 0.68; the resulting adjust- 
ments in solar units are M^/L^ = 0.76 M*/Ly, and 
M*/Lij = 1.33 M*/L/. For the maximum-mass dark 
matter halo, our dynamical measurement of M^/L^ = 
4.6^ ; 7 M /L 0i _r agrees with population-based measure- 
ments from lCappellari et al.l (|2006l ). who find M*/L/ « 3 
M /L 0! / (M*/L fi ft 4 M /L „r) for M87 and several 
other early-type galaxies observed with SAURON. How- 
ever, a different range of stellar mass-to-light ratios i s 
suggested by the study of Ivon der Linden et al.l (|2007[) . 
who model the stellar masses of 625 BCGs from SDSS. 
The peak values of their mass and luminosity distribu- 
tion functions yield a ratio of M*/Ly = 3.1 M^/L^.y 
MJL R = 2.4M /L , fl ). Additionally, [G raves k Fabcr 
2010) have compiled various estimates of M*/L for a 
sample of ~ 16, 000 early-type galaxies from SDSS, find- 
ing that M*/L increases with <r, and M*/Ly ft 2 — 3 
M Q /L 0! v (M*/L R ft 1.5-2.3 M /L 0! _r) for the highest- 
dispersion objects in their sample (a ~ 250 km s" 1 ). 
Both of these results fall significantly below our dynam- 
ical values of M*/Lfl for NGC 6086. The discrepancy 
is most severe for orbit models without dark matter 
(M*/L fi = 6.7±° J M /L 0)B ). 

To test whether any model of the gravitational poten- 
tial for NGC 6086 can hold lower values of M+/Lr and 
still fit our data, we have run a series of orbit mod- 
els with M+/Lr fixed at 2.5 M /L 0j ij;. We sampled 
the LOG dark matter halo parameters v c and r c over 
a wide range of values and marginalized y 2 over trials 
with M, = 3.5 x 10 9 M and M, = 7.0 x 10 9 M . The 
resulting best-fit parameter values are v c = 420 km s _1 
and r c = 2.0 kpc. Fixing v c and r c at these values and 
finely sampling M, yields M. = 5.5±o g x 10 9 M , includ- 
ing systematic errors, with Xmin = 1018.6. In compar- 
ison, Xmin — 1010.2 for our maximum-mass LOG halo 
(v c = 500 km s _1 ; r c = 8.0 kpc). We conclude that 
the lower values of M*/Lr motivated by stellar popula- 
tion modeling can produce a reasonable fit to our kine- 
matics, given a more centrally concentrated dark matter 
halo. Nonetheless, our original assumptions about the 



dark matter halo profile of NGC 6086 produce a better 
fit than the assumption of low M*/L^. Decreasing the 
enclosed stellar mass leads to a larger best-fit black hole 
mass, matching the trend from our other trials. 

For each individual dark matter halo, the dominant 
systematic effect in our results is a discrepancy in the 
wings of LOSVDs from IFS versus long-slit data, illus- 
trated by a jump in from negative to positive values. 
The magnitude of the resulting error is 22% in M, and 
9% in M^/Lfl. Uncertainties in determining an optimal 
template star and AO PSF yield smaller errors, totaling 
~ 9% in M, and ~ 2% in M*/Lfl. Additional, unmea- 
sured systematic errors may arise from our assumptions 
of axisymmetry and edge-on inclination. 

Our investigation of NGC 6086 has established that 
stellar-dynamical measurements of M, in BCGs are pos- 
sible with existing facilities; this dramatically expands 
the sample volume of viable targets. We plan to follow 
this work with similar measurements from an ongoing 
survey of several additional BCGs. This survey is a 
critical step toward a statistically robust census of black 
holes in the Universe's most massive galaxies, and will 
eventually shed new light on the histories of galaxies 
and black holes at the hearts of galaxy clusters. 
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APPENDIX 

A: UNCERTAINTIES IN THE PSF AND STELLAR TEMPLATES FOR OSIRIS 

We ran additional series of orbit models to assess the effects of stellar template mismatch and PSF uncertainty 
on the best-fit values of M. and rvL/Lft. We expect these systematic effects to be largest for OSIRIS data. Lower 
signal-to-noise spectra have larger uncertainties in template matching, and the crowded series of atomic and molecular 
features in ii-band is more sensitive to template choice than the Call features at 0.85 /im. The structure of the 
LGS-AO PSF is sensitive to atmospheric turbulence, laser power, the density and thickness of the ionospheric sodium 
layer, performance of the wavefront sensor and deformable mirror, and the brightness and position of the tip/tilt star. 
Because several of these factors change over time, delays between observations of an extended science target and a point 
source can induce errors the estimated PSF. Uncertainties in seeing-limited PSFs arise predominantly from changes in 
atmospheric turbulence. In order to highlight the effect of template and PSF errors for OSIRIS, we excluded the more 
stable GMOS data from our trials below. Each trial included our maximum mass dark matter halo. We summarize 
the trials in Table [3] 

In Figure [El we display various estimates of the average PSF for OSIRIS observations of NGC 6086. Our initial 
estimate, displayed at the top middle, was constructed from a one-time sequence of exposures of the LGS-AO tip/tilt 
star. We have characterized the two-component structure of this PSF by fitting a narrow Gaussian profile plus a 
broad Moffat profile. The resulting FHWM values are 0.10" and 0.42". The narrow component contains 44% of the 
total flux. However, we cannot directly compute a Strehl ratio from this percentage, because the 0.05" pixels in our 
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TABLE 3 
PSF and Template Star Trials 



Data 


PSF 


Template 


v c 


(UpC) 


M. 

V iU m Q) 


M+/L R 

v iVi © / ^© J 


x- 2 ■ 


Ndof 


(1) 




(2) 


(3) 


(4) 


(5) 


(6) 


(7) 


(8) 


(9) 


OSIRIS + 


CBH99 


A 


M4III 


500 


8.0 




4 -°-0.6 


419.0 


570 


OSIRIS + 


CBH99 


B 


M4III 


500 


8.0 


l.9I Ll 


4 6+°' 7 


412.7 


570 


OSIRIS + 


CBH99 


B 


MOIII 


500 


8.0 


2 0+ 13 


4 5+0.7 


415.6 


570 


OSIRIS + 


CBH99 


C 


M4III 


500 


8.0 


9 1+ 1 - 4 


4 -°-0.6 


414.6 


570 



Notes: Column 1: Data sets included in trial. Column 2: estimated PSF for OSIRIS data. A: original PSF measured from 
tip/tilt star with OSIRIS, folded over major and minor axes of NGC 6086 (Figure [151 top right). B: tapered and folded PSF 
(Figure [151 bottom left). C: core-halo PSF, with 25% Strehl ratio (Figure [T5l bottom right). Column 3: spectral type of 
template star for OSIRIS data. We used spectra from HD 110964 (M4III) and HD 108629 (MOIII). Column 4: circular velocity 
of LOG dark matter halo (Eq. [3]). Column 5: core radius of dark matter halo (Eq. [3]). Column 6: best-fit black hole mass. 

Column 7: best-fit i?-band stellar mass-to-light ratio. Column 8: minimum \ 2 value for all models. Column 9: degrees of 
freedom in model fits to LOSVDs. Computed values include a smoothing factor of 1 degree of freedom per 2 velocity bins for 
non-parametric LOSVDs from OSIRIS. Quoted errors in M, and IVR/L^ correspond to 68% confidence intervals. All trials in 
Section H used PSF B and template star HD 110964 (M4III). 

collapsed OSIRIS images undersample the .ff-band diffraction limit (X/D = 0.033"). Our trials use PSFs derived from 
the original tip/tilt images, rather than the Gaussian plus Moffat model. The PSFs used in our trials (A, B, and C) 
were all symmetrized with respect to the major- and minor-axis position angles of NGC 6086, in order to match the 
spatial folding of kinematic data for axisymmctric models. 

PSF A (Figure [151 top right) symmetrizes the estimate from the tip/tilt star, but contains no further changes. This 
estimate has two primary limitations. It only captures the PSF at a single moment in time, and it does not account for 
the 42" separation between the tip/tilt star and the center of NGC 6086. A complementary method for estimating the 
PSF is to compare the collapsed OSIRIS mosaic of the galaxy center to an image from HST /WFPC2. We convolved 
this mosaic with an appropriate WFPC2 PSF and convolved the WFPC2 image with a trial PSF, so that both images 
would have the same total smoothing kernel. Using our original PSF estimate from the tip/tilt star, we found that 
the brightness profile of the convolved WFPC2 image was shallower than the collapsed OSIRIS mosaic. Therefore, we 
tested a series of PSFs that multiplied this original estimate by tapering functions with various inner and outer radii. 
We found the lowest root-mean-squared residual for PSF B, depicted in Figure [TBI (bottom left), which suppresses 
power at radii beyond 0.3". In order to probe the effect of a Nyquist-sampled, diffraction-limited core, we constructed 
a third estimate, PSF C (Figure [TBI bottom right). We interpolated the original tip/tilt-based PSF from 0.05"/pixel 
to 0.01"/pixel, and replaced the central 0.1" x 0.1" with a 2-dimensional Airy profile (X/D — 0.033"). The amplitude 
of the core was scaled relative to the outer profile to yield a Strehl ratio of 25%. 

PSFs A, B, and C incorporate a variety of methods for PSF estimation, and the range of variation in their structures 
is comparable to the uncertainty in any particular method. We find a standard d eviation of 8% in M, and 1% in 
M*/L/{ among trials using each of these three PSFs. Similarly, iNowak et al.l (|2008h found little variation in M. and 
M+/Lk s from different estimates of the PSF for AO data of Fornax A. For all of our primary trials to meas ure M, 
and M+/L R in NGC 6086 (SectionH]), we matched OSIRIS data to PSF B. For data from GMOSand [CBH99I, we use 
Gaussian PSFs. The CMOS PSF has a FWHM of 0.4", based on acquisition camera images. [CBHH report a FWHM 
of 1.3". Uncertainties in these seeing-limited PSFs should have smaller effects on measurements of M,, due to poorer 
spatial resolution. 

Our stellar template library for OSIRIS contains stars with spectral types from G8 through M4. In principle, we can 
determine a best-fit weighted combination of library templates for each spectrum of NGC 6086, while simultaneously 
fitting for the LOSVD. However, noise in our OSIRIS spectra renders this method unstable: when we include a 
large number of templates, the weights vary dramatically between adjacent spatial regions. To narrow our sample, 
we compared the equivalent widths (EWs) of jff-band absorption features in NGC 6086 to each library template, as 
illustrated in Figure [16] EWs for NGC 6086 were measured from a high-S/N spectrum covering a large portion of 
the OSIRIS field-of-view. Our kinematic fitting is dominated by high-EW features at 1.54 /jm, 1.55 /jm, 1.62 /im, and 
1.64 /im. For these features, the three M-giant stars provide much closer EW matches than other stars in our library. 
Trial kinematic fits using only M-giant stars weighted HD 110964 (M4III) by nearly 100%; for simplicity, we adopted 
this star as our only template for the final kinematic extraction from OSIRIS spectra. 

To test template mismatch, we extracted a second set of LOSVDs from OSIRIS spectra, this time using HD 108629 
(MOIII) as the template star. The difference between spectral types MOIII and M4III is a reasonable estimate for 
template mismatch in our kinematic fitting, as indicated by our equivalent width measurements. In particular, our 
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Fig. 15. — PSFs estimated for OSIRIS observations of NGC 6086. Top left: horizontal slice through the center of each PSF. Top middle: 
Original estimate from the LGS-AO tip/tilt star. Top right: PSF A, folded for axisymmetry. Bottom left: PSF B, folded and tapered. 
Bottom right: PSF C, interpolated for a diffraction-limited core and 25% Strehl ratio. The major and minor axes of NGC 6086 are oriented 
45° from the plotted X- and Y-axes. Relative flux is defined with respect to the peak value of PSF C. All of our orbit models used PSF B 
with OSIRIS data, with the exception of trials to test model dependence on PSF structure. 

fits exclude the Mg/Fe absorption feature near 1.50 rest, which is a better match to our K-dwarf template. Orbit 
model trials with M4III-based LOSVDs versus MOIII-based LOSVDs differ by 5% in M. and 2% in M*/L R . 

Even though HD 110964 compares favorably to other template stars, its average equivalent width over the spectral 
features used to extract LOSVDs is ~ 13% higher than the average equivalent width of NGC 6086. To improve the 
EW agreement, we artificially decreased the line strengths of the template spectrum, by a constant factor Jew = 0.83. 
We did not allow Jew to vary over different spatial regions in NGC 6086. Optimizing the value of Jew produces 
smoother LOSVDs, which are also more stable to small changes in spectral smoothing. However, we have found a 
positive correlation between Jew and the value of derived from the resulting fit. In non-parametric terms, the wings 
of the best-fit LOSVD become truncated as absorption features in the template star are artificially mad e shallow er. 
This effect may contribute partially to the differences between our computed values of /14 and those from[CBH99. 
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Fig. 16. — Equivalent widths of NGC 6086 and template stars for different H-band spectral features. The shaded areas mark the spectral 
range used for kinematic extraction from OSIRIS spectra. Template stars are: HD 110964 (M4III); HD 108629 (M0III 1); HD 63348 (M0III 
2); HD 44537 (K5I); 15 Lib (K4III); HD 99492 (K2V); HD 94386 (K2III); 55 Cnc (G8V); HD 89638 (G8III). 
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B: NON-PARAMETRIC LOSVDS FROM OSIRIS AND CMOS 



We present our final ex tracted LOSVDs from OSIRIS (Figure TB) and CMOS (Figure M) , and LOSVDs derived 
from the measurements of lCBH99l (Figure HSJ . The LOSVDs from OSIRIS use template star HD 110964 (M4III). We 
compare each LOSVD to the best-fitting orbit model with our maximum-mass LOG halo, and the best-fitting orbit 
model with no dark matter. Including dark matter gives a better fit to the total set of LOSVDs, with a cumulative 
difference A\ 2 = 104.2. Similarities between the LOSVDs generated by each model occur in part because ~ 30, 000 
orbits give the models a high degree of flexibility to optimally fit the data. The models inability to perfectly match 
the data arises in part from constraints such as axisymmetry and uniformity in M*/Li?. 
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Fig. 17. — LOSVDs from OSIRIS spectra. Solid blue lines are corresponding LOSVDs generated by the best-fitting orbit model with 
the maximum-mass LOG dark matter halo (M+fLn = 4.7 Mq/Lq, M. = 3.5 X 10 9 Mq, v c = 500 kpc, r c = 8.0 kpc). Dotted red lines 
are generated by the best-fitting orbit model with no dark halo (Mi/L^; = 6.8 Mq/Lq, M, = 3.2 X 10 8 Mq). In each sub-plot, r is the 
distance from the center of NGC 6086, and 8 P is the range of angles with respect to the major axis. Negative values of r indicate spectra 
from the south side of NGC 6086. For each LOSVD, A% 2 is the difference in the \ 2 statistic for the two models: A\ 2 > indicates that 
the model including dark matter (solid blue line) is a better fit. 
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Fig. 18. — LOSVDs from GMOS spectra. Solid blue lines are corresponding LOSVDs generated by the best-fitting orbit model with 
the maximum-mass LOG dark matter halo (M^/Lfl = 4.7 M Q /L , M. = 3.5 X 10 9 M , v c = 500 kpc, r c = 8.0 kpc). Dotted red lines 
are generated by the best-fitting orbit model with no dark halo (M*/L/j = 6.8 Mq/Lq, M, = 3.2 X 10 s Mq). In each sub-plot, r is the 
distance from the center of NGC 6086, and 6 P is the range of angles with respect to the major axis. Negative values of r indicate spectra 
from the south side of NGC 6086. For each LOSVD, A\ 2 is the difference in the x 2 statistic for the two models: Ax 2 > indicates that 
the model including dark matter (solid blue line) is a better fit. 
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Fig. 19. — LOSVDs generated from the kinematic moments measured by CBH99. Solid blue lines are corresponding LOSVDs generated 
by the best-fitting orbit model with the maximum-mass LOG dark matter halo (M*/L fl = 4.7 M /L Q , M. = 3.5 X 10 9 M , v c = 500 kpc, 
r c = 8.0 kpc). Dotted red lines are generated by the best-fitting orbit model with no dark halo (M*/L^j = 6.8 Mq/Lq, M. = 3.2 X 10 8 
Mq). In each sub-plot, r is the distance from the center of NGC 6086, along the major axis. Negative values of r indicate spectra from 
the south side of NGC 6086. For each LOSVD, Ax 2 is the difference in the \ 2 statistic for the two models: A\ 2 > indicates that the 
model including dark matter (solid blue line) is a better fit. 
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